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Abstract 

A nonlinear internal model control (NIMC) strategy that 
incorporates the nonlinear model structure and the estimator 
dynamics in the control law is presented for the control of 
complex dynamical systems characterized by input-output 
multiplicities, nonlinear oscillations and chaos. A model 
based estimator is designed to provide the unmeasured 
process states that capture the fast changing nonlinear 
dynamics of the process to incorporate in the controller. The 
estimator uses the mathematical model of the process in 
conjunction with the known process measurements to 
estimate the states. The design and implementation of the 
estimator supported NIMC strategy is studied by choosing 
two typical continuous non-isothermal nonlinear processes, 
a chemical reactor and a polymerization reactor, which show 
rich dynamical behaviour ranging from stable situations to 
chaos. The results evaluated under different conditions show 
the superior performance of the estimator based NIMC 
strategy over the conventional controllers for the control of 
complex nonlinear processes. 
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Introduction 

Control of nonlinear systems exhibiting complex 
dynamic behaviour is a challenging task because such 
systems present a variety of behavioural patterns 
depending on the values of their physical parameters 
and intrinsic features. Among these systems, 
continuous chemical reactors and continuous 

polymerization reactors have received special 

attention due to their complicated dynamics and 
economic importance. Depending on the parameter 
values, these systems can be operated at steady state 
or present oscillatory and chaotic motions. Several 
researchers have analysed continuous stirred tank 
chemical reactors to characterize the phenomena of 
multiple steady states, simple and modulated 


oscillations and chaos [Uppal, Ray and Poore, 1976, 
Schmitzet al., 1979, Doedel and Heinemann, 1983, 
Blanco and Bandoni, 2007]. Certain continuous 
polymerization reactors exhibit highly nonlinear 
dynamics because of the complicated reaction 
mechanisms associated with the large number of 
interactive reactions. Many researchers have 
investigated the existence of steady state multiplicities; 
self-sustained oscillations and chaos in continuous 
solution polymerization reactors [Jaisinghani and Ray, 
1977, Schmidtet al., 1984, Teymour and Ray, 1989, 
Flores-Tlacuahuac et al., 2005]. 

The nonlinear dynamical systems due to variations in 
the operating parameter space can lead to exhibit 
open-loop unstable behaviour and some chemical 
systems are more sensitive to these variations where 
small perturbation in the operating conditions can 
cause the stable operation to degenerate into 
oscillatory and chaotic motion. The unstable, 
oscillatory and chaotic phenomenon displayed by the 
chemically reacting systems has desirable as well as 
undesirable features. The desirable feature of multiple 
steady state condition is that one of the unstable 
steady state may correspond to higher rate of 
reaction/yield and this becomes the preferred state to 
be operated on it for enhancing the process 
performance. The desirable feature of chaos is that it 
enhances mixing and chemical reactions and provides 
a vibrant mechanism for transport of heat and mass. 
On the other hand, the intrinsic features of the reacting 
systems with the interactive influence of chemical or 
thermal energy may cause irregular dynamic 
behaviour leading to degraded performance. In such 
situations, chaotic behaviour is considered as 
undesirable and should be avoided. Chaotic processes 
show extreme sensitivity to initial conditions and the 
process trajectories can diverge exponentially making 
the task of control more difficult. However, chaos 
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offers great flexibility to operate chemical systems 
because there are an infinite number of unstable 
periodic orbits (UPO) embedded in a chaotic attractor, 
in which one can choose a specific UPO along with its 
time averaged pre-specified performance defined in 
terms of conversion, yield or selectivity and the 
process can be stabilized to operate it on the chosen 
UPO. Various conventional and advanced techniques 
have been reported for the control of open-loop 
unstable chemical and polymerization reactors and the 
reactors that show oscillatory and chaotic dynamics. 
As far as unstable reactors are concerned, Choi and 
Ray [1988] and Chan et al. [1993] have employed 
standard PID controllers to stabilize open-loop 
unstable steady state conditions of an olefin 
polymerization reactor, however, Takamatsu et al. 
[1986] and Kwag and Choi [1994] showed that such 
controllers may also pose serious control problems if 
the polymerization reactor presents significant 
nonlinear behaviour. Dittamar et al. [1991] reported 
that linear predictive controllers might lead to proper 
control of open-loop unstable polymerization reactor. 
Various nonlinear model predictive control strategies 
have been successfully applied to the control of open- 
loop unstable CSTR systems [Hidalgo and Brosilow, 
1990, Hapoglu, 2002, Venkateswarlu and Venkat Rao, 
2005]. Nonlinear processes characterized with chaotic 
dynamic behaviour can be controlled either by 
stabilizing the UPOs in the neighbourhood of the 
desired UPO or stabilizing the dynamic trajectory of 
the system exactly at the unstable steady state. Various 
controllers based on ad-hoc, conventional and 
advanced approaches have been reported to stabilize 
the chaotic systems. In the ad-hoc approach, the 
properties of chaos are used to develop ad-hoc 
methods without considering conventional control 
techniques. The most commonly used method is the 
OGY method by which one of the unstable orbits in 
the chaotic attractor can be stabilized via small time 
dependent perturbations to a system parameter [Ottet 
al., 1990, So and Ott, 1995]. The potential difficulty 
associated with the OGY method is the long transient 
times before the system enters local region where the 
perturbation is effective [Petrov et al., 1994]. In 
conventional control approach, various methods 
including a proportional-integral (PI) controller 
[Pellegrini and Biardi, 1990] and a modified PI/PID 
controller [Bandyopadhyay, et al., 1997] have been 
employed for controlling and operating the chaotic 
reactors under favourable conditions. However, the 


complex nature of nonlinear dynamical systems 
severely limits the use of conventional linear 
controllers to provide the desired operating 
performance. In advanced control approach; different 
nonlinear model predictive controllers have been 
reported to stabilize the oscillatory and chaotic 
dynamics in continuous chemical and polymerization 
reactors [Lima, 1996, Qammaret al., 1996]. 

Nonlinear control is a class of advanced control in 
which the nonlinear process model serves as the basis 
for controller design. This type of control is expected 
to provide improved control performance since the 
control structure preserves the nonlinearities of the 
real process in the form of its mathematical model. The 
generic model control (GMC) introduced by Lee and 
Sullivan [1988], the globally linearizing control (GLC) 
proposed by Kravaris and Chung [1987] and the 
nonlinear internal model control (NIMC) proposed by 
Henson and Seborg [1991] are the prominent nonlinear 
model based controllers. The GMC allows the 
implementation of the nonlinear process model 
directly into the controller structure. The GLC is a 
model based controller that transforms a nonlinear 
input/output system into a linear input/output system 
through a nonlinear transformation. The NIMC 
proposed by Henson and Seborg is different from the 
GMC of Lee and Sullivan [1988] and the GLC of 
Kravaris and Chung [1987] in that it includes implicit 
integral action in the control structure by using the 
difference between the process output and model 
output as a feedback signal. A nonlinear filter is 
employed in NIMC which provides a tuning 
parameter that can be adjusted for process/model 
mismatch. Since reasonably accurate mathematical 
models have been developed to predict the 
complicated dynamics of nonlinear systems, design of 
nonlinear controllers based on these models is a useful 
control alternative for complex dynamical systems. 
This work is focused towards the model based control 
of nonlinear dynamical systems in the frame work of 
NIMC approach. The NIMC approach based 
controllers have been reported earlier for nonlinear 
systems where the dynamic behaviour is not so 
complex [Kurtz and Henson, 1997, Venkateswarlu and 
Gangiah, 1997]. 

The objective of this paper is to derive a controller 
based on NIMC approach for complex dynamical 
systems that incorporates the nonlinear model 
structure and the estimator dynamics in the controller 
formulation. A model based estimator is designed to 
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provide the unmeasured process states that capture 
the fast changing nonlinear dynamics of the process to 
incorporate in the controller. The estimator uses the 
mathematical model of the process in conjunction with 
the known process measurements to estimate the 
states. The design and implementation of the 
proposed estimator based NIMC strategy is studied by 
choosing two typical continuous chemical and 
polymerization reactors that present challenging 
operational and control problems due to their complex 
open-loop dynamics such as input-output 
multiplicities, parametric sensitivity, nonlinear 
oscillations and chaos. Further the comparison of the 
estimator based NIMC strategy is made with 
conventional controllers. 

Control Algorithm 

The general form of a single input single output (SISO) 
system with state space description is: 


e = y d -(y-y) ( 6 ) 

where y is the process model output. The feedback 
signal is simplified to e=yd with perfect model 
assumption. A nonlinear filter is employed in NIMC 
which provides a tuning parameter that can be 
adjusted for process/model mismatch. The new input v 
for NIMC is defined as 

v = -TL r ;h(x)-T r __ ] L'~ 1 h(x)-.... 

- v x h(x) + z x e 

where v are controller tuning parameters; r is the 
relative order and e=yd. According to this approach, 
the control law for the system is given by 

V - L f r h(x) 

U = — = p(x) + q(x)v (8) 

L g L f h(x ) 

Application Processes 


X = fix) + g(x)u 
y = h{x) 


Two typical continuous nonlinear dynamical 
processes are chosen for the design and 
(2) implementation of the estimator based NIMC strategy. 


where x is the vector of states, u is the manipulated 
input, y is the measured output, f(x) and g(x) are 
vector functions and h(x) is the scalar function. The 
relative order of the system defined by Eqs. (1) and (2) 
is expressed by the following relations: 

LgL f k h{x) — 0, k < r — \ (3) 

L g L f h(x) * 0 ^ 

where r represents the relative order of the system and 
Lfh(x) is the Lie derivative of the scalar function h(x) 
with respect to the vector function f(x) with L°f 
h(x)=h(x). Similarly higher order Lie derivatives as well 
as the Lie derivative of the scalar function h(x) with 
respect to the vector function f(x), and then with 
respect to the vector function g(x) can be defined. The 
relative order defined by Eqs. (3) and (4) represents 
the number of times the output y must be 
differentiated with respect to time so that the input u 
appears explicitly. 

A general form of the control law is written as 
u=p(x) +q(x) v (5) 

where v is the new input. 

The NIMC approach includes an implicit integral 
action by using the difference between the plant and 
model outputs as a feedback signal: 


Che mica l re actor 


A non-isothermal, irreversible, first order series 
reaction A— »B— >C in a CSTR with control input and 
load disturbances is described by the following 
dimensionless mass and energy balance equations 
[Kahlertet al., 1981]: 


dx i 

= 1 — Xj - Da x i exp(x^ j\ + £ j x 2 ) + dy 


dx^ 

- = -^2 + Da x i exp(x^ j\ + s 

- Da S.x^expfCr^ / 1 + £^ 4 X 3 ) + 

= -jcj +B Da x 1 exp(x 2 j\+e ^* 3 ) 

- Da B aSx 2 exp(kx^ j\ + £jX^ ) - P(x^ - ) (11) 

+ /? u f + d ^ 


( 9 ) 

( 10 ) 


The variables xi and xi denote the dimensionless 
concentrations of species A and B respectively, and X3 
is the dimensionless temperature. The parameter X3c 
represents the reactor coolant temperature. An 
externally manipulated variable lit can be defined to 
denote a measure of the deviation in the coolant 
temperature from the reference value X3c. The load 
disturbances in feed compositions are denoted by di 
and cfe, and the load disturbance in the reactor 
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temperature is denoted by cfo. This reactor system 
exhibit multi-stationary behaviour, oscillations and 
chaos for the parameter values shown in Table 1. 

TABLE 1 CHARACTERIZATION OFSTEADYSTATES 


Solvent Mass Balance 



P S f V sf 


Ps 


+P S V S 


dT d(l/p s ) 
dr dT 


(14) 


Set 

No 

Parameter values 

XI s 

X2s 

X3s 

Stability 

i 

D a=0 .0 6,S=0. 000 5, eA=0, 
k=l, a=0.426, (3=7.7,B=55.0 

0.0378 

0.9501 

6.0500 

Unstable 

n 

Da=0.26,S=0.5,sA=0,k=l, 
a=0.426, p=7.7,B=57.77 

0.0729 

0.1259 

3.8902 

Stable 

limit 

cycle 

HI 

Same as Set Ilexcept 
P=7.9999 

0.0819 

0.1391 

3.7627 

Chao tic 


Homopo lymerization reactor 

Continuous homopolymerization reactors represent 
one important class of dynamical systems. The 
homopolymerization reactions are characterized by 
high heat release, high viscosity and poor heat transfer. 
Vinyl acetate homopolymerization has shown 
complex dynamic behaviour for certain ranges of 
operating conditions due to the large heat of 
polymerization, the gel effect and the large activation 
energy of the initiation step. These factors can readily 
produce phenomena such as multiple steady states, 
steady state instability, limit cycles and chaos. 


Initiator Mass Balance 


dc ; 


, L = c if -cMoVi + k d°) 
dr 

Energy Balance 

dT P / (T f - T) 0 ah 


(15) 


+ 


dr p pC p (MW) n 

0 UA(T-T C ) 


V mPm k n P 


(16) 


pc> 


V 


+ P u , 


hr the above equations, q Q q i is a term that takes into 
account changes in the density of the reactive mixture 
and is defined by 


qSi = — + — + mW) mPm R„ 

Pm Ps 


Pm 

\Pp ) 


(17) 


<TT_ 

dr 


( 


Ps V s 


dQ „ d (\ //O „ d(\jp p ) 


dT 


+ P m v, 


dT 


+ Pp v p 


dT 


and /3 — 


9 UA 

~pC p ~V~ 


The overall density and specific heat of the mixture are 


The kinetic mechanism for free-radical vinyl acetate 
homopolymerization reaction is as follows. 

/ — — — » 2 R (Radical Initiation) 


R + M- 

~^P l 

(12) 

P t +M- 

K " > P 
' 1 1 + 1 

(Chain Propagation) 

P,+Pn~ 

— 

i+n 

(Chain Termination) 


The mathematical model for homopolymerization 
reaction of vinyl acetate in tertiary butanol using 
AIBN as initiator is described by the following 
equations [Pinto and Ray, 1995]: 

Monomer Mass Balance 


dv m _ 

dr 

+ V ,nP 


P mf^mf 

— — -~q 0 9i 

Pm 

dT d(l/p a ) 


(MW) m R m 9 

Pm 


(13) 


P = Pm V m +Ps V s +Pp V p 

c p =y m Cpn,+y s c P s+y P c PP 

Under quasi steady state assumption for the free- 
radical species (R and P), the normalized rate of 
monomer reaction is expressed as 

Rate m = k p P + R ckr (18) 

where concentration of species P can be written as 

p = R ‘ k ^ M (19) 


The R k is the normalized rate of initiation as given 
by 


R 


ckr 


M 


( 20 ) 


By substituting the value of Rckr from Eq. (20) in Eq. 
(18), we have 


4 


dr dT 




AdvancedChemical Engineering Research (ACER) Volume 2Issue 1, March 2013 


ww w.s e ipub .org/ace r 


Rate m =k p P + 




M 

RJMW) m 


P V 

m n 


(21) 


where R is the rate of monomer consumption 
expressed by 


R. 


k„Pp m v m 

{MW) m 


+ 2 fk d Ci 


( 22 ) 


The termination and propagation constants ( k t ,k p ) 

used in the above equations include gel effect with the 
form 



(23) 


The gel effect correlations are 
Sp =1 

g, = exp(-0.4407x, - 6.753x, 2 -0.3495*, 3 ) 


In Eq. (23), k ; ° and k p are the kinetic constants at zero 
polymer concentration, g t andg^ are the gel effect 
TABLE 2 VA HOMOPOLYMERIZATION SYSTEM PARAMETERS 


f= 0.8 

V = 500 (ml) 

UA = 12 (cal/min K) 

AH = 2 100 0(c al/g mo 1) 

Cpm = 0.470 (cal/g K) 

Cps = 0.716 (cal/g K) 

p m (T) = 958.4 - 1 .3276(7’ - 273) (g/lt) 


p(T) = 1211- 0.8496(T - 273) (g/lt) 


P,{T) = 


74120 


(60.21 + 0.1167) 


(g/lt) 


C pp =0.3453 + 9.55 xlfT 4 (T -298) (caV g K) 

-6100 

k° p = 82.212 xl0 8 e^ (it/ g moi.s) 

-30800 


k d = 94.8 x 10 L ^ e RT (it/gmoi.s) 

-2462 

k* =469.392x10 10 (lt/gmoLs) 


correlations that take into account the effects of 
increasing polymer concentration. The other 
parameters necessary for simulation are presented in 
Table 2. 


Design of Estimator Based Nonlinear 
Controller 


The nonlinear internal model controller (NIMC) relies 
on mathematical model of the process incorporating 
process state variable information in the controller 
formulation. Since the state variables desired by the 
controller are not easily available through 
measurement or available with large measurement 
delays, a model based nonlinear estimator is designed 
to provide the unmeasured states needed by the 
controller. 

Controller Des ign 

The design of NIMC for complex dynamic processes 
considered in this work is briefed as follows: 


1) Chemical Reactor 


The f(x) and g(x) in Eq. (1) can be written by 
omitting the load disturbances from the 
dimensionless mass and energy balances in Eqs. (9) 
- ( 11 ): 


fix)- 


1 - Xj - Dax i exp| 


l + £ A x 3j 


- *2 + DaXy exp| 


~ x n+ BDax^ exp| 


1 + £ A x 3 


- DaS*2 ex P| 


kx. 


1 + £ , x . 
A 3 


f 1 


f kx, ) 

3 

- BDaaSXr. exp 

3 

1 + £ 

1 + £ 

2. 1 

^ A 3) 


'y a 3 J 


-px 


g(x) 


0 

0 


p 


(24) 

(25) 


The controlled output in Eq. (2) is given by 

y - h{x) = X 3 (26) 

From Eqs. (3) and (4), the relative order of the 
system is evaluated as one. Thus for r=l, the new 
input v in Eq. (7) is obtained as 

V = r, (*3 - X 3 ) (27) 

On computing the Lie derivatives and substituting 
them along with the new input in Eq. ( 8 ) leads to 
the following control law 
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r, (x 3 — x 3 ) — 


— x 3 + BDax l exp 


u . = - 


x 3 

V 1 + s a x 3 j 


— BDaaSx 2 exp 


kx 3 

yl + e A x 3 J 


- px 3 


P 


where n is the NIMC tuning parameter. 
2) Polymerization reactor 


(28) 


Th ef(x) and g(x) in Eq. (1) can be defined from the 
mass and energy balances of Eqs. (13) - (16): 


f(x) = 


W Wm {MW) ~ R 


p e e 

m 

P sf V sf Ws 


m m | y p f 'm 1 


m m Jr 


p fVf T 


9 AH 9 UA(T-T) 

v p k P — 

pC (MW) m m p pC y 


+ p v 


m / p s ) 


p 6 9 s s dT 


if % q i 


p A T f- T 

p 


9 AH , n 9 UA(T-T) 

+ v p k P — 

pC (MW) m m p pC y 


— - k,c. 

9 9 1 d 1 


P fl T f~ T ) K Hv m P m k p P UA(T~T c ) 


9p 


pC (MW) V pC 
p in p 


g(*) = 


(30) 

The controlled output in Eq. (2) is obtained as 
y = h(x) - T (31) 

The conditions in Eqs. (3) and (4) defining the 


(29) 


polymer reactor system is of relative order one. The 
new input v in Eq. (7) is computed as 


v = r 1 (J f -r) 


(32) 


f 


t i( T s -T)- 


Pf (Tf ~ T) AHv m k pPm P UA (T-T) 


dp pC p (MW) m V pC 

J3 


p J 


(33) 


Estimator Design 

The success of estimator based controller for nonlinear 
dynamical systems relies on the performance of the 
estimator as well as the model of the process that 
supports the estimator. Since mathematical models 
with enough details are used to develop the NIMC 
strategy for complex dynamical systems, it is 
advantageous to use the same models in the estimator 
design in order to obtain fast and accurate estimation 
of unmeasured process states desired by the controller. 
Recently, model based estimators have been 
successfully employed for state estimation in various 
stable systems [Venkateswarlu and Gangiah, 1992, 
Schuler and Schmidt, 1993, Sargantanis and Karim, 


1994, Aguilar-Lopez and Martinez-Guerra, 2005]. hr 
this work, a nonlinear model based estimator known 
as extended Kalman filter (EKF) is designed to obtain 
the state estimates from the known temperature 
measurements of the respective processes. The 
estimator design for NIMC strategy of this work is 
similar to that of the recently reported estimator 
involved in the GLC control of a chaotic chemical 
reactor [Karri et al., 2009, Karri, 2011]. The general 
process representation for model based state 
estimation is given in Appendix A. The EKF 
estimation algorithm is given in Appendix B. The 
elements of the state transition and measurement 
matrices involved in the estimator are given in 
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Appendix C . 

Results and Discussion 

Chemical Reactor 

The reactor system described by Eqs. (9)-(ll) exhibits 
multi-stationary behaviour, oscillations and chaos for 
parameter values given as sets I, II and III, respectively, 
in Table. 1. For parameter set I, the system shows 
unstable steady state, where as for parameter sets II 
and III, the system exhibits limit cycle oscillations and 
chaotic behaviour. The temperature data used for state 
estimation are obtained through numerical integration 
of model equations using Gear's method with a 
sampling time of 0.0001 units. A model based state 
estimator, extended Kalman filter is employed to 
estimate the states xi, X 2 and X3 using the temperature 
measurement of the reactor. The performance of the 
EKF estimator is evaluated by considering different 
cases of process parameter values shown in Table 1. 
The results in Fig. 1 represent the phase plane plots of 
the actual and estimated states of CSTR corresponding 
to the parameter values in set III. 


applied with the objectives of controlling the chemical 
reactor at an unstable steady state in the multiplicity 
region, at the unique unstable steady state responsible 
for limit cycle oscillations and at unique unstable 
steady state responsible for chaotic motion. The 
controller parameter involved in NIMC is tuned as 
n=14.0 based on objective function to minimize the 
integral of absolute error (IAE). The results are also 
compared with those of a modified PID controller 
with feedback mechanism, which has the form given 

by [22], 

~ s of (0 e + s of + s i\ (34) 

dt dt J 

where f(t)=t and/'(f)=l. The tuning parameters of the 
PID controller (ffl,a) are evaluated using Ziegler - 
Nichols method. A sampling time of 1 sec is used for 
the implementation of the estimator and the controller. 
For parameter values in set I, the system exhibits 
multiple steady state behaviour. The controller goal is 
to shift the process operating at an arbitrary point (xio 
= 0.04, X 2 o = 0.9, X3o =5.75) to the set condition in Table 1 
and maintain it at that state. 




FIG. 1. ACTUAL (Xj , X 2 ) AND ESTIMATED (xj , X 2 
PHASE PLANE PLOTS FOR SET - HI 


Fig. 2 shows the actual and estimated concentration 
profiles corresponding to set El in the presence of 
stochastic noise in temperature obeying Gaussian 
distribution. Similar results are evaluated for set I and 
set II. In all the cases, the estimated states are in close 
resemblance with the actual states. 

The NIMC strategy supported by the state estimator is 


0.16 -I 

0.14 - 
0.12 - 
0.1 -I 
* 0.08 - 
0.06 - 
0.04 - 
0.02 
0 -- 

23456789 10 

Time [mill] 


0.16 



0 -I ■ , , . , , . , 

23456789 10 


Tune [mill] 

FIG. 2ACTUAL (xj ) AND ESTIMATED (x : ) 
CONCENTRATION PROFILES FOR SET-HI IN THE PRESENCE 
OF STOCHASTIC DBTURBANCE 

Fig. 3 shows the process output and controller output 
plots of NIMC and PID controller, respectively for this 
condition. For parameters in set II, the system exhibits 
sustained oscillatory behaviour (stable limit cycle) and 
the controllers are also required to regulate the 
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trajectory at the unique unstable steady state, which is 
the desired condition. Fig. 4 shows the results of the 
NIMC and PID controller to realize the set condition in 
Table 1 from an arbitrarily selected initial state (xio = 
0.08, X20 = 0.103, X3 o= 3.654). For parameters in set III, the 
system exhibits chaotic behaviour. In this objective, the 
controller has to stabilize the chaotic trajectory exactly 
at the corresponding unique unstable steady state. To 
realize this objective, the controllers are also applied 
using the same initial condition as in set II to satisfy 
the desired condition for set III given in Table 1. 

The process output and controller output plots of 
NIMC and modified PID controller corresponding to 
this objective are shown in Fig. 5. These plots indicate 
that the controllers fulfil the desired objectives without 
any offset. The nonlinear controller has shown better 
stabilization over PID controller. The performance of 
the controllers is also evaluated by applying them for 
controlling the system at the desired conditions in the 
presence of deterministic and stochastic load 


disturbances. A deterministic disturbance of 1.0 and a 
stochastic load disturbance generated through random 
Gaussian noise of zero mean and a standard deviation 
of 0.25 are considered to represent the d3 in 
temperature measurement. The results of the 
controllers for stabilizing the system at unstable steady 
state responsible for chaotic motion with these 
disturbances conditions are shown in Figs. 6 and 7. 

These results indicate the better performance of the 
NIMC strategy in the presence of either type of load 
disturbances. The results of the controllers are also 
evaluated for a series of step changes in the set points 
of the controlled variables. The results in Fig. 8 show 
the process output and controller output plots of 
NIMC and modified PID controller for parameters in 
set III for a series of ±20% step changes in the set 
conditions in Table 1. These results show the better 
performance of the estimated supported nonlinear 
controller. 


r 


5 


0 


(a) NIMC 


2 3 4 5 

Time [min] 



Time [min] 


NIMC 



PID 



Time [min] 


FIG. 3 PROCESS AND CONTROLLER OUTPUT PLOTS FORSET I OF (A) NIMC AND (B) PID 
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NIMC 


PID 



FIG. 8 PROCESS AND CONTROLLER OUTPUT PLOTS OF NIMC AND PID FORA SERIES OFSTEP CHANGES IN SET POINTS OF THE 

CONTROLLED VARIABLES FORSET m 


Polymerization Reactor 

The mathematical model of homopolymerization 
reactor is solved by numerical integration using Gear's 
method with a sampling time of 0.0001 units. The 
bifurcation diagram of VA homopolymerization 
shown in Fig. 9 drawn for the operating conditions, Tc 
= 318 K, Tf= 315 K, vmf= 0.3 and c//= 0.03203 gmol/lt, 
presents the existence of multiple steady states, limit 
cycle and bifurcation points as well as stable and 
unstable zones with respect to the variation of 
residence time. 

From Fig. 9, it can be observed that the transition 
occurs from stable to unstable zone beginning at 0 = 24 
min. This transition causes oscillatory behaviour in 
response that eventually lead to period doubling 
bifurcation cascade and chaos as shown in Fig. 10(a). 
The magnification results of Fig. 10(a) are further 
explained in Figs. 10(b) and 10(c), which show the 
sequence of periodic windows of periods 2, 4, 8 and so 
on leading to chaotic region. At 0 = 27.37 min, the 
chaotic behaviour with huge temperature oscillations 
is observed as shown in Fig. 11. The time domain 
response is shown in Fig. 11(a) and the corresponding 


phase plane plot is depicted in Fig. 11(b). 



FIG. 9 PERIOD DOUBLING BIFURCATION DIAGRAM 

A model based estimator, extended Kalman filter (EKF) 
is employed to estimate the states Vm vs a and T using 
the temperature measurement of the reactor. The 
performance of the estimator is evaluated under 
different conditions. Fig. 12 compares the actual and 
estimated results of VA homopolymerization reactor. 
These results confirm the usefulness of the method of 
EKF as a state estimator for polymerization reactor. 
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Residence time(min) 


FIG. 10 ONE DIMENSIONAL POINCARE BIFURCATION 
DIAGRAM: (A) TRANSITION FROM LIMIT CYCLES TO CHAOS, 
(B) MINIMUM TEMPERATURE AND (C) MAXIMUM 
TEMPERATURE 
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FIG. 11 CHAOTIC BEHAVIOUR (A) TIME DOMAIN RESPONSE, 
(B) PHASE PLANE PLOT 






Temperature (K) Temperature |K) 

FIG. 12 ACTUAL AND ESTIMATED PROFILES OF MONOMER VOLUME FRACTION (A) TIME DOMAIN RESPONSE (B) PHASE PLANE 

PLOTS 
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FIG. 13 PROCESS AND CONTROLLER OUTPUT PLOTS OF NMC AND PI 


The NIMC strategy supported by the state estimator is 
applied with the objectives of controlling the 
homopolymerization reactor under different 
conditions and the results are presented in terms of 
three cases, hr the first case, the process is initially 
operated under stable operation. Then changes are 
introduced in temperature set point and the controller 
is applied to maintain the desired operation under 
different input disturbance conditions introduced at 
different timings. In the second case, the process is 
operated under oscillatory condition. The controller is 
applied to stabilize the process at the desired 
condition at a time during the open-loop oscillatory 
dynamics, hr the third case, the process is operated 
under chaotic behaviour. The controller is applied to 
suppress the chaotic dynamics and to maintain the 
desired operation. The results of the NIMC strategy 
are also compared with those of a conventional PI 
controller. The tuning parameter in NIMC strategy is 
chosen as n =16.5. The PI controller parameters are 
evaluated by using Ziegler and Nichols method and 
further tuned and set as Kc = 1.76 and n =4.0. A 
sampling time of 1 sec is used for implementing the 
estimator and controller. For the first case, the process 
is initially at the stable steady state with T= 366.4 K and 
Tc= 320 K. At time t = 8.3 min, the temperature set-point 
is changed from 366.4 to 364.87 K and at t = 200 min, Tf 
is changed from 315.00 to 313.00 K while keeping the 
set-point as 364.87 K. The input and output profiles 


corresponding these conditions are shown in Figs. 13(a) 
and 13(b), respectively. In the same case, at t = 400 min, 
the monomer volume fraction, Vmf is changed from 0.3 
to 0.33 while keeping the set-point at 364.87 K and Tf at 
313 K. The input and output profiles corresponding 
these conditions are shown in Figs. 13(c) and 13(d), 
respectively. Though the changes are introduced at 
different timings in one time scale, the results are 
represented in separate figures with two time scales 
for the sake of clarity. From these results, it is observed 
that the NIMC tracks the desired condition very 
quickly, whereas the PI controller exhibits oscillatory 
behaviour for a longer duration. The PI controller 
could even fail for change in Vmf as shown in Fig. 13(c). 
The NIMC effectively rejects the disturbance 
conditions and leads the process back to its desired 
operation in a very conservative fashion. 

hr the second case, the process is disturbed from its 
initial steady state by changing the Tc from 320 K to 
318 K at f=8.3 min. This condition makes the process to 
exhibit oscillatory behaviour. This open-loop 
behaviour of the system is shown in Fig. 14(a). During 
the oscillatory operation at f = 217.4 min, the 
controllers are applied to stabilize the operation and to 
maintain the desired response condition, which is set 
as 364.87K for this case. The results of process output 
and controller output for both the controllers are 
shown in Figs. 14(b) and 14(c), respectively. From 
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these results, it is observed that both the controllers 
are able to track the process to the desired condition, 
but the NIMC has shown better performance by 
quickly suppressing the oscillations. In the third case, 
the process is exhibiting chaotic behaviour for T=347.6 
K and Tc = 318 K. 




FIG. 14 OPEN LOOP AND CLOSED LOOP RESULTS: (A) OPEN- 
LOOP BEHAVDUR (B) CONTROLLEDOUTPUT(C) CONTROL 
INPUT 

During the chaotic operation at t = 240 min, the 
controllers are applied to stabilize the operation and to 
maintain the desired response which is set as 364.87K 
for this case. The results of process output and 
controller output for both the controllers are shown in 
Figs. 15(a) and 15(b), respectively. Again, the NIMC 
has performed better by quickly suppressing the 
chaotic dynamics with less stringent control actions. 

The results evaluated for both the chemical and 
polymerization reactors show the superior 
performance of the estimator based NIMC strategy 
over the conventional controllers. The performances of 
the controllers for different conditions are also 
expressed in terms of IAE values. The IAE values 


shown in Table 3 correspond to the conditions for 
which the results are depicted in figures referred in 
Table. 3. These quantification results further confirm 
that the NIMC strategy is better suited for the control 
of nonlinear dynamical processes. 




FIG. 15 OPEN LOOP AND CLOSED LOOP RESULTS: (A) 
CONTROLLED OUTPUT (B) CONTROL INPUT 


TABLE 3 CONTROLLER'S PERFORMANCE 


Chemical Reactor 

IAE 

REFERENCE 

NIMC 

PID 

0.938 

365.163 

Fig. 3 

0.748 

72.763 

Fig. 4 

0.347 

93.524 

Fig. 5 

7.121 

83.042 

Fig. 6 

0.924 

93.709 

Fig. 7 

6.852 

118.350 

Fig. 8 

Polymer Reactor 

0.000985 

174.193 

Fig. 13 

40.605 

110.697 

Fig. 14 

138.626 

339.425 

Fig. 15 


Conclusions 

Complex dynamical systems can present the 
phenomena such as multiple steady state behaviour. 
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nonlinear oscillations and chaos. This phenomenon is 
usually analysed by means of dynamic mathematical 
models. The availability of such models can form a 
basis for the development of advanced nonlinear 
controllers for such systems. The nonlinear model 
based controllers require the process state information 
in the controller formulation. This necessitates the 
need for an estimator to provide the unmeasured 
process states desired by the controller. Thus, this 
work is focused towards the development of a 
nonlinear model based controller for complex 
dynamic systems, which is supported by a nonlinear 
model based estimator. A nonlinear internal model 
control (NIMC) strategy that incorporates the 
nonlinear model structure and the estimator dynamics 
in the control law is presented for the control of 
complex dynamical systems. The estimator is designed 
to provide the unmeasured process states that capture 
the fast changing nonlinear dynamics of the process to 
incorporate in the controller. The performance of the 
estimator supported NIMC strategy is evaluated by 
applying it for the control of a non-isothermal 
nonlinear chemical reactor and a homopolymerization 
reactor, which exhibit rich dynamical behaviour 
ranging from stable situations to chaos. The results 
evaluated under different conditions show the 
superior performance of the estimator based NIMC 
strategy over the conventional controllers. 

Appendix A: Process Representation 

The time varying model of the nonlinear process is 
represented by 

x(t) = + w(t) , x(0) = X 0 (A.l) 

where x(t) is n dimensional state vector, / is a 
nonlinear function of state x(t) and w(t) is an additive 
Gaussian noise with zero mean. The nonlinear 
measurement model with observation noise can be 
expressed as 

y(t k ) = h(x(t k )) + v(t k ) (A.2) 

where h is a nonlinear function of state x(tk). The 
expected values of noise covariance matrices for the 
initial state x(0), process noise w(t) and observation 
noise v(tk) are given by the following relations, 

P 0 = £[(-*o -*(0))O 0 - V(0)) r ] 

0(O = £[w(t)w T (t)] (A.3) 

R(t k ) = EM t k )v T (t k )] 


where Pois initial statecovariance matrix, Q (t) is 
process noise covariance matrix and R(fo) is 
observation noise covariance matrix. The EKF 
estimation algorithm is given in Appendix B. 


Appendix B: State Estimation Algorithm 

The extended Kalman filter is computed in two steps. 
The first is a prediction step, which is used to 
extrapolate the previous best estimates, and the 
second is a correction step by which the updated 
estimates are formed. Since prediction is based on 
process model, continuous prediction and discrete 
correction are employed in the estimation scheme. 

Prediction equations: By starting with an initial estimate 
Xo and its covariance P o at time zero and no 
measurements are taken between tk-i and tk , the 
propagating expression for the state estimate and its 
covariance from tk- i to tk are. 


x(t/t k -i) = f(x(t/t kA ),t) (B.l) 

Ht/t k - 1 ) = F (x(t/t kA ), t )P(t/t kA ) 2 

+ P{t/t k _ x )F T (x{t/t k _ x ),t) + Q{t) 


wher eF(x(t / t k _ x ),t) is the state transition matrix 
whose ij th element is given by 


F(x(t/t k _ x ), t) 


dXj(t ) 


x(t)=x(t/t kA ) 


(B.3) 


The solution of the propagated estimate x(t/t k _ x ) and 
its covariance P(f/fc-i) at time fc are denoted by 
x{t k /t k .\) and P(tk/tk- i). By using measurements at time 

tk, the update estimate x(t k /t k ) and its covariance 
P(tk/tk) are computed. 

Correction equations: The equations to obtain corrected 
estimates are: 


x(t k /t k ) = x(t k /t k _ l ) + K(t k )[y(t k )-h(x(t k /f t _j))] 

(B.4) 

P(t k /t k ) = (I-K(t k )H(x(t k )))P(t k /t k _ x ) (B.5) 


K(t k ) = 


P(t k /t k _ x )H T (x(t k )) 


(H(x(t k ))P(t k /t k _ x )H T (x(t k )) + R) 


(B.6) 


where, H(x(t k )) = 


dhj(x(t k )) 

dx(t k ) 


At k )=x(t k /t k _ i) 


The recursive initial conditions for state and 
covariance are: 
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x(t k /t k _ i) = x(t k /t k ) 


Appendix C: Elements ofF and Ft Matrices Chemical 
reactor 

The nonlinear dynamic model of the reactor system in 
its dimensionless form is used in conjunction with the 
temperature measurements to estimate the reactor 
species concentrations. The design of soft sensor 
involves the following components: 

State vector 

The chemical species concentrations and the 
temperature define the state vector as 

x = [X\X 2 x }t ] T (C.l) 

State transition matrix 

The elements fij of the state transition matrix. Fare 
computed by taking the partial derivatives of f(x) 
defined by (24) - (26) with respect to the state vector: 


fn 

fn 

fn 

h 1 

f 22 

f 2 3 

/31 

fn 

fn 


(C.2) 

Measurement matrix 

The measurement relation for temperature is 

H =[0 0 1] (C.3) 

The temperature state equation, (26) in its discrete 
form is used as the nonlinear measurement equation, 
h(x). The elements of the measurement transition 
matrix. Hx are computed by taking the partial 
derivatives of h(x) with respect to the state vector 

H x — [h\ ] h l2 ^ 13 ] (C-4) 

All these components are evaluated for chaotic reactor 
and used with the EKF estimator to obtain measured 
and unmeasured states of the reactor. The soft sensor 
uses the temperature data of every sampling time as 
measurements and provides the estimates of 
temperature as well as reactor species concentrations. 

Polymer reactor 

The nonlinear dynamic model of the VA 
homopolymerization system is used in conjunction 
with the temperature measurements to estimate the 
volume fraction of monomer, solvent, concentration of 


initiator and reactor temperature. The design of soft 
sensor involves the following components: 

State vector 

The volume fraction of monomer, solvent, 
concentration of initiator and reactor temperature 
define the state vector as 

* = [ V ,n v , C i T ( C5 ) 

State transition matrix 

The elements fij of the state transition matrix. Fare 
computed by taking the partial derivatives of f(x) 
defined by (28) - (30) with respect to the state vector: 


fn 

fn 

fn 

fu 

fn 

fll 

fn 

f 2A 

fn 

fn 

fn 

fu 

u 

fn 

fn 

fu 


Measurement matrix 

The measurement relation for temperature is 

H = [0 0 0 1] (C.7) 

The temperature state equation, (30) in its discrete 
form is used as the nonlinear measurement equation, 
h(x). The elements of the measurement transition 
matrix, Hx are computed by taking the partial 
derivatives of h(x) with respect to the state vector 

H x —[^11 ^12 ^13 ^14] (C-8) 

All these components are evaluated for 
homopolymerization reactor and used with the EKF 
estimator to obtain measured and unmeasured states 
of the reactor. The state estimator uses the temperature 
data of every sampling time as measurements and 
provides the estimates of temperature as well as 
reactor species concentrations. 

Notation 

B heat of reaction parameter 

a initiator concentration in reactor, gmol/l 

a/ initiator feed concentration, gmol/l 

Cp heat capacity of reaction mixture, cal/g °C 

Cpf heat capacity of feed stream, cal/g °C 
di, d2, d.i load disturbance in xi, X 2 , X3 
Da Damkohler number 

e set point error 
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/ initiator dissociation efficiency 

-A Hr heat of polymerization, cal/gmol 

k ratio of the activation energies for the series 

reaction 

kd rate constant for initiator decomposition. l/s 
kp rate constant for propagation, l/gniol s 
kt rate constant for termination, l/gmol s 

Kc controller gain 

(MW)m molecular weight of monomer, g/gmol 
P live radical concentration, gmol/l 

q 0 qi ratio of outlet to inlet volumetric flow rate 
Rm rate of consumption of monomer, gmol/l min 

5 ratio of the rate constants for the series 
reaction 

t time, min 

T reactor temperature, K 

Tc coolant temperature, K 

Tf feed temperature, K 

UA overall heat transfer coefficient, cal/ min K 

in manipulated variable 

V volume of the reactor, l 

xt mole fraction of polymer 

xi dimensionless concentration of species A 

X 2 dimensionless concentration of species B 

X3 dimensionless temperature 

X3c dimensionless coolant temperature 

j 2 set , xa^set point for variable xi, X 2 and X3 

xis , X 2 s , xassteady state values of xi, X 2 and X3 

y; weight fraction of respective components 

(i=m,s,p) 

Greek symbols 

a ratio of heat effects for the series reaction 
P heat transfer coefficient 

so initial controller gain 

si initial controller gain 

sa activation number 

6 residence time, min 


z dimensionless time (t/0) 

id time constant for derivative action 

v time constant for integral action 

Vm volume fraction of monomer in reactor 

Vmf volume fraction of monomer in feed 

p density of reaction mixture, g/l 

pf density of feed, g/l 

pur) density of feed adjusted to reactor 

temperature, g/l 

p m density of monomer at reactor conditions, g/l 
p m f density of monomer at feed conditions, g/l 

Subscripts 
m monomer 

s solvent 

p polymer 

/ feed 
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